  #include "ms_spongebob.h"
#include "ms_voxelimage.h"
#include "ms_cylinder.h"
#include "ms_misc.h"
#include <magick/MagickCore.h>
#include <string.h>
#include <stdlib.h>
#include <stdio.h>
#include <math.h>

//--- #include "unistd.h"
#include "ms_distribution.h"

#include "ms_distribution.cxx"
#include "ms_sphere.cxx"
//--- #include "fe_vecmat.cxx"
#include "ms_misc.cxx"
#include "ms_voxelimage.cxx"



#ifndef PI
#define UNDEF_PI
#define PI 3.1415926535897932384626433
#endif

using namespace std;

/* ---------------------------------------------------- */

Range SpongeBob::validateRandomValues(RandomDistribution* useDistribution, const int count) 
{ 
  
/*        int nOkay = 0;
        int nLowerZero = 0; 
        int nHigherOne = 0;
*/	
        double random;
	
	
	Range range;
	range.lowest=1.0; range.highest=0.;
	
        for(int i = 0; i < count; i++)
        { 
	  random = useDistribution->GetRandomValue();
          if(range.lowest > random) range.lowest = random;
          if(range.highest < random) range.highest = random;
/*          if(random < 0.0)
            nLowerZero++;
          if(random > 1.0)
            nHigherOne++;
          else
            nOkay++;
*/	  
        } // next i
	
/*	cout << "Tested " << count << " random values" << endl;
       cout << nLowerZero << " values < 0.0; " 
            << nOkay << " values okay; " 
            << nHigherOne << " values > 1.0; " << endl;
       cout << "lowest value = " << fLowest 
            << "; heighest value = " << fHighest << endl;*/
//	range.lowest=fLowest; range.highest=fHighest;
	return range;
	

} // END OF validateRandomValues()



/* ********************************************** */
// Test für die charakteristische Zahl
// Einfacher 2x2x2 oder 3x3x3 Würfel wird befüllt und charakteristische Zahl ausgegeben
// Vergleich mit von Hand ausgerechneten Zahlen
// z.B.Würfel mit 4 befüllten voxeln: char. Zahl: 1020


/* int main(int argc, char *argv[]){
  
        VoxelImage vx_img( 1, 3, '\0' );
	char fullname[1024];
	
        vx_img.SetData(0,0,0);
        vx_img.SetData(1,0,0);
        vx_img.SetData(0,0,2);
        vx_img.SetData(2,0,2);
//        vx_img.SetData(1,1,0);
//        vx_img.SetData(1,0,1);
//        vx_img.SetData(0,1,1);

	vx_img.ShowNumber(fullname);

	return 0; 
  
 
}  */


 int main(int argc, char *argv[]){

    //---if (argc!=3)
	//---{
    //---    return 0;
	//---}

	bool bPartialOverlapping = false;
	bool bSameCores = false ;
	bool bConstRadius = false;//+++

	double cv[13] = { // Parameter der Verteilung
	0.3 , 0.75 , 1 , 1.25 , 1.5 , 1.75 , 2 , 2.5 , 2.75 , 3 , 
	4 , 6 ,  10	
	};
	 

	
	for (int l=2 ; l<13; l++)
	{
	for (int k=3 ; k<13; k++)
	{
		double Parameter1 = cv[l];
		double Parameter2 = cv[k];
		double minRadius = 0.005;
		double maxRadius = 0.035;
		for(int i=0;i<20 ;i++) // wieviele Schritte, wie voll bepacken?
		{ // start i
		  int S_bin = 1000 + (500*i);
		  double edgeLength = 1.;

		cout << "Anzahl der Kugeln: " << S_bin << endl;  
		SpongeBob Sbob( S_bin , edgeLength );


		if((argc > 1) && (strcmp(argv[1], "true" )==0)) bPartialOverlapping = true;
		if((argc > 2) && (strcmp(argv[2], "true" )==0)) bSameCores = true;
		if((argc > 3) && (strcmp(argv[3], "true" )==0)) bConstRadius = true; 

  // 	minRadius = strtod(argv[3],NULL);
  // 	maxRadius = strtod(argv[4],NULL);
		Sbob.setVerbose(false); // Ausgaben an cout machen
		//---RandomDistribution* RadiusDistribution = new ExponentialDistribution( mu ); 
		RandomDistribution* RadiusDistribution = new BetaDistribution( Parameter1, Parameter2 ); //+++		
		sprintf(Sbob.CurrentDistributionName , "BetaDistribution");
		sprintf(Sbob.CurrentParameterName, "%f-%f", cv[l] , cv[k]);
		
		RandomDistribution* PlaceDistribution = new UniformDistribution(-0.5,0.5);
		PlaceDistribution->setSeed(k*l);
		RandomDistribution* CoreDistribution = new BetaDistribution( Parameter1, Parameter2 );
		
		Sbob.createNewSpheres( RadiusDistribution, PlaceDistribution , CoreDistribution , S_bin, minRadius, maxRadius, bPartialOverlapping, bSameCores, bConstRadius);
		
		
		delete RadiusDistribution; // Distribution wieder loeschen
		delete PlaceDistribution;
		delete CoreDistribution;

		if(Sbob.CancelSbin == true) {
  		  cout << "Parameter1 " << cv[l] << endl;
		  cout << "Parameter2 " << cv[k] << endl;

		  cout << "abgebrochen bei Kugel " << S_bin << endl;
		  break;
		}
		
		  
		} // end i
	} // end k
	} // end l

	return 0; 
}  

SpongeBob::SpongeBob( const int n_spheres, const double edgeLength ){
/* Konstruktor */
	S_bin	= n_spheres; // im Header definfiert
	pSpheres = new sphere [n_spheres];
		// erstellen der Kugeln im Raum
		// createNonOverlappingSpheres(RandomDistribution* RadiusDistribution, int S_bin, double minRadius, double maxRadius);
		// VoxelImage vx_img( 1, 200, ' ' );
	fEdgeLength = edgeLength;
	pDistributionUsed = NULL; //
	pDistributionRadius = NULL; //
	CoordDistributionUsed = NULL; //
	minRadius = 0.01; //
	maxRadius = 0.02; //
	bVerbose = false; // do no output to cout
	CancelSbin = false; //
}

SpongeBob::~SpongeBob()
{
	delete [] pSpheres;
	pSpheres	= NULL;
	S_bin		= 0;
	fEdgeLength	= 1.0;
}

////////////////////////////////////////////////////////////////////////////////
//+++ Die Funktion createRadiusSpheres(...) unterscheidet sich von createSpheres() 
//+++ nur darin, dass in einem fall ein konstanter Radius fuer alle Spheres benutzt 
//+++ wird. Das laesst sich leicht ueber einen Parameter bConstRadius steuern.

/** 
  * Diese Funktion erzeugt die Poren in dem SpongeBob.
  * @Param RadiusDistribution enthaelt einen Zeiger auf die zu benutzende Verteilung
  *        fuer den/die Radien der Poren
  * @Param S_bin, die Anzahl der Poren. Wenn S_bin ausserhalb kleiner gleich 0 oder
  *        groesser dem S_bin Wert des SpongeBob ist, werden die maximal moeglichen 
  *        Poren fuer den SpongeBob erzeugt.
  * @Param minRadius, der Minimalwert fuer den Radius der Poren.
  * @Param maxRadius, der Maximalwert fuer den Radius der Poren.
  * @Param bPartialOverlapping ist true wenn sich die Poren teilweise ueberlappen duerfen.
  * @Param bSameCores true, wenn die Kerne der Poren gleich gross sind
  * @Param bConstRadius ist true, wenn alle Poren mit dem selben Radius erzeugt werden.
  * @Return true, wenn alle Poren unter den angegebenen Parametern erzeugt werden konnten.
 */
bool SpongeBob::createNewSpheres(RandomDistribution* RadiusDistribution, RandomDistribution* PlaceDistribution , RandomDistribution* CoreDistribution , 
								 int S_bin, 
								 double minRadius, double maxRadius, bool bPartialOverlapping, 
								 bool bSameCores, bool bConstRadius) 
{	bool bSuccess = true; // Rueckgabewert
	// Parameter Uebernehmen:
	this->pDistributionUsed = RadiusDistribution;
	this->CoordDistributionUsed = PlaceDistribution;
	this->CoreDistributionUsed = CoreDistribution;
	this->minRadius = minRadius;
	this->maxRadius = maxRadius;
	
	char fullname[1024]; //fn fuer imagestack
	char fullname2[1024]; //fn fuer vtk
	char fullname3[1024]; //fn fuer vtk
	
	
	if((S_bin <= 0)  && (S_bin > this->S_bin)) //*** MAX_SPHERES ersetzt durch this->S_bin
		S_bin = this->S_bin;				   //*** MAX_SPHERES ersetzt durch this->S_bin
	
	const int nMaxTry = 30; // maximale Versuche eine Kugel i zu erzeugen
	int nTry = 0;
	int Sfail = 0;
	
	double fCoordinates[3] = {0.0,0.0,0.0}; 
	double fRadius = 0.;
	double fRadiusCoreMin = 0.005; // Minimum, den der Kern haben darf
	double fRadiusCore = 0.;//--- = RadiusDistribution->GetRandomValue() * (fRadius - RadiusCoreMin) + RadiusCoreMin; // Kern der Pore.
	
	
	Range rangeradius = validateRandomValues(pDistributionUsed); //+++ hier werden fHighest und fLowest ausgegeben
	Range rangecore = validateRandomValues(CoreDistributionUsed);
	

	if(bConstRadius)
	{	fRadius = 0.03;
	
	// getRandomRadius(1, minRadius, maxRadius, rangeradius); // fRadius = getRandomRadius();
		if(fRadius <= 0.0) // Suche abgebrochen
			bSuccess = false;
		else if((bPartialOverlapping == false) && (fRadius >= fEdgeLength)) // hier kann nur gleichheit auftreten
			bSuccess = false; // aber beruehren ist im Fall von nicht ueberlappend verboten.
	}
	
	if(bSuccess && bSameCores)
	{	fRadiusCore = getRandomRadius(1000, fRadiusCoreMin, fRadius, rangecore);
		if(fRadiusCore <= 0.0) // Suche abgebrochen
			bSuccess = false;
	}
	const int nCloneCount = 27;
	sphere periodicSpheres[nCloneCount]; //*** s_out -> periodicSpheres

	VoxelImage vx_img( 1, 200, '\0' );

	for(int i = 0; ((i < S_bin) && (!CancelSbin)); i++)
	{	nTry = 0; // Soll verhindern, dass bei einer schlechten Verteilung ewig versucht wird neue Poren zu erstellen.
		do // {...} while((bSuccess == false) && (nTry < nMaxTry)); 
		{	bSuccess = true; // until other is proofen
			
			
			CoordDistributionUsed->GetRandomVector(fCoordinates, 3); // neuer Mittelpunkt
			//if(range.highest > 1. || range.lowest < 0.)  // Werte aus validateRandomValues
			//{	
			//	if(bVerbose)	cout << "Coordinates (x; y; z;): ";
			//	for(int k = 0; k < 3; k++)
			//	{ // normieren der Koordinaten auf den Wuerfel und versetzen in die Mitte des Wuerfels
			//	
			//		fCoordinates[k] = (fCoordinates[k] - range.lowest)*( 1 / ( range.highest - range.lowest ) ); // verschieben der Koordinaten auf die Länge von 0-1
			//		fCoordinates[k] = (fCoordinates[k] - 0.5) * fEdgeLength;
			//		if(bVerbose)	cout << fCoordinates[k] << "; ";
			//	} // next k
			//	if(bVerbose)	cout << endl;
			//}
			//else
			//{
			//	if(bVerbose)	cout << "Coordinates (x; y; z;): ";
			//	for(int k = 0; k < 3; k++)
			//	{ // normieren der Koordinaten auf den Wuerfel und versetzen in die Mitte des Wuerfels
			//		fCoordinates[k] = (fCoordinates[k] - 0.5) * fEdgeLength;
			//		if(bVerbose)	cout << fCoordinates[k] << "; ";
			//	} // next k
			//	if(bVerbose)	cout << endl;
			//}

			
			if(bConstRadius == false)
			{	fRadius = getRandomRadius( 1 , minRadius , maxRadius, rangeradius ); // nur ein Versuch fuer den Radius, weil sowieso 1000 mal versucht wird, Pore erneut zu erstellen.
				if(fRadius <= 0.0) // Suche abgebrochen
					bSuccess = false;
				else if((bPartialOverlapping == false) && (fRadius >= fEdgeLength))
					bSuccess = false;
			}

			if(bSuccess && bPartialOverlapping)
			{ 	if(bSameCores) // testen ob Kern im Radius liegt, falls (bConstRadius == false)
				{	if(fRadiusCore > fRadius)
						bSuccess = false;
				}
				else
				{	fRadiusCore = getRandomRadius( 1 , fRadiusCoreMin , fRadius , rangecore);
					if(fRadiusCore <= 0.0) // Suche abgebrochen
						bSuccess = false;
				}
			}


			if(bSuccess) // Radius (und ggf. Kern) ist okay
			{	pSpheres[i] = sphere(fCoordinates, fRadius); // Pore erstellen
				pSpheres[i].coreR = fRadiusCore;
				if(bVerbose)	cout << "created sphere index i = " << i << "; try = " << nTry << "; R = " << pSpheres[i].R << "; coreR = "<< pSpheres[i].coreR << endl;
				PeriodicCopy( pSpheres[i] , periodicSpheres ); 	// Kopien erstellen
				// Bedingungen fuer jede periodische Kopie der neuen Pore i gegen alte Poren j testen:
				for(int nCloneID = 0; ((nCloneID < nCloneCount) && bSuccess); nCloneID++)
				{	// jeden Klon gegen alle Poren j testen.
					//*** Fehlt: Noch Besonderheit bei Radien >= fEdgeLength beachten!!!
					//***		Wenn (bPartialOverlapping == true),
					//***        dann muss j ==i auch getestet werden, ausser nCloneID == 4;
					//***        Bei (bPartialOverlapping == false) kollidieren sie,
					//***         wird aber schon beim Radius erstellen abgefangen.
					for (int j = 0; ((j < i) && bSuccess); j++) // alle bisherigen Poren durchgehen
					{	// Die folgenden if(...)-Abfragen koennte man zu einem einzigen if() zusammenfassen.
						// Ein einziges if() waere aber noch schwerer nachzuvollziehen.
						if(bPartialOverlapping == true) // ueberlappend
						{	double fMidDistance = PointDistance(periodicSpheres[nCloneID].X, pSpheres[j].X);
							if(		(fMidDistance < (periodicSpheres[nCloneID].R + pSpheres[j].coreR) ) // Kern wird von Pore ueberlagert
								||	(fMidDistance < (pSpheres[j].R + periodicSpheres[nCloneID].coreR) ) // beidseitig 
							  ) // Ende Bedingungen
							{	nTry ++; // wegen Kollision, weiterer Versuch eine Pore i zu erstellen.
								bSuccess = false; // Schleifen j und nCloneID abbrechen
							} // end if
							//--- else if(bVerbose)	cout << "Success at j = " << j << "; clone " << nCloneID << endl;
						} // END OF if(bPartialOverlapping == true)
						else // nicht ueberlappend
						{	if(periodicSpheres[nCloneID].Collision(pSpheres[j])) // Klon kollidiert mit alter Pore
							{	nTry ++; // wegen Kollision, weiterer Versuch eine Pore i zu erstellen.
								bSuccess = false; // Schleifen j und nCloneID abbrechen
							} // end if
							//--- else if(bVerbose)	cout << "Success at j = " << j << "; clone " << nCloneID << endl;
						} //end else bPartialOverlapping
						
						if((bSuccess == false) && bVerbose)
						{	cout << "Fail at sphere j = " << j << " with clone " << nCloneID << " (x; y; z;) = ";
							for(int k = 0; k < 3; k++)
								cout << periodicSpheres[nCloneID].X[k] << "; ";
							cout << " r = " << periodicSpheres[nCloneID].R << endl;
						}
					} // NEXT j
				} // NEXT nCloneID
			} // END if(bSuccess) // Radius ist okay
			if(bSuccess)
				{ for(int nCloneID = 0; ((nCloneID < nCloneCount) && bSuccess); nCloneID++)    
					periodicSpheres[nCloneID].DrawVoxelImage(vx_img, 'Z');
				}
// 				pSpheres[i].DrawVoxelImage(vx_img, 'Z');

		} while((bSuccess == false) && (nTry < nMaxTry)); // end do
	if(nTry >= nMaxTry) Sfail++;
	if (Sfail >= 1) CancelSbin = true;
	// cout << "i: " << i << endl;
	} // NEXT i


//	cout << "Sfail: " << Sfail << endl;
      
	sprintf(fullname,"%f/", fRadiusCoreMin);
	if(mkdir( fullname , 0777)==-1) cerr<<"Error: " <<strerror(errno)<< endl;
	
	sprintf(fullname,"%f/%s/", fRadiusCoreMin, CurrentDistributionName);
	if(mkdir( fullname , 0777)==-1) cerr<<"Error: " <<strerror(errno)<< endl;

	sprintf(fullname,"%f/%s/%f/", fRadiusCoreMin, CurrentDistributionName , maxRadius);
	if(mkdir( fullname , 0777)==-1) cerr<<"Error: " <<strerror(errno)<< endl;

	sprintf(fullname,"%f/%s/%f/%f/", fRadiusCoreMin, CurrentDistributionName , maxRadius , minRadius);
	if(mkdir( fullname , 0777)==-1) cerr<<"Error: " <<strerror(errno)<< endl;

	sprintf(fullname,"%f/%s/%f/%f/%i/", fRadiusCoreMin , CurrentDistributionName , maxRadius , minRadius , S_bin);
	if(mkdir( fullname , 0777)==-1) cerr<<"Error: " <<strerror(errno)<< endl;
	
	sprintf(fullname,"%f/%s/%f/%f/%i/%s", fRadiusCoreMin , CurrentDistributionName , maxRadius , minRadius , S_bin, CurrentParameterName);
	if(mkdir( fullname , 0777)==-1) cerr<<"Error: " <<strerror(errno)<< endl;


	strcpy(fullname2, fullname);
	strcpy(fullname3, fullname);
	
	strcat(fullname2, "/spheres.vtk");
	strcat(fullname3, "/Stack");

	vx_img.WriteVTKFile(fullname2);
	vx_img.WriteImageStack(fullname3);
	
	
	
	sprintf(fullname , "%f/%s/%f/%f" , fRadiusCoreMin , CurrentDistributionName , maxRadius , minRadius );
	strcat(fullname , "/Data.csv");

	// Einfügen der Daten in eine .csv Datei
	ofstream scriptfile;
	scriptfile.open( fullname , ios::app|ios::out );
	scriptfile << CurrentDistributionName << ";" << S_bin << ";" << maxRadius << ";" << minRadius << ";" << CurrentParameterName << ";";
	
	scriptfile.close();

	
	vx_img.ShowNumber(fullname);


 	return bSuccess;

} // END OF bool SpongeBob::createNewSpheres(RandomDistribution* RadiusDistribution, int S_bin, double minRadius, double maxRadius, bool bPartialOverlapping, bool bSameCores, bool bConstRadius) 

////////////////////////////////////////////////////////////////////////////////
/** Radius fuer Pore (oder Porenkern) ermitteln.
 * @Param min kleinster Wert fuer Radius
 * @Param max groesster Wert fuer Radius. Wenn max <= 0.0 werden min und max aus den 
 *        Attributen minRadius und maxRadius genommen.
 * @Return gefundener Radius, oder Wert kleiner gleich 0, wenn mit der gewaehlten
 *         Verteilung kein Radius gefunden werden konnte.
 */
double SpongeBob::getRandomRadius(const int nMaxTrys, double min, double max , Range range)  //+++ wegen Verteilungen mit Wertebereichen ausserhalb von [0,1]
{ 	double fRadius = 0.0;
	if(max <= 0.0)
	{	min = minRadius;
		max = maxRadius;
	}
	RandomDistribution* pDistribution;
	if(pDistributionRadius != NULL) // besondere Verteilung fuer Radien?
		pDistribution = pDistributionRadius;
	else
		pDistribution = pDistributionUsed;
	if(pDistribution != NULL) // irgendeine Verteilung ist initialisiert
	{	
//		Range range = validateRandomValues(pDistribution); //+++ hier werden fHighest und fLowest ausgegeben
		int nTry = 0;
		bool bRadiusOkay = false;
		do // {...} while( (bRadiusOkay == false) && (nTry++ < nMaxTrys));
		{	
			if(range.highest > 1. || range.lowest < 0.)  // Werte aus validateRandomValues
			{	
				fRadius = pDistribution->GetRandomValue(); // Zufallswert rauspicken
				fRadius = (fRadius - range.lowest)*(1/( range.highest - range.lowest ) ); // ihn auf die Strecke von 0-1 trimmen		
				fRadius = fRadius * (max - min) + min;
				if(	  (fRadius > 0.0) // Kugel mit keinem oder negativem Radius, ist keine Kugel! Auf wenn minRadius was anderes sagt.
				   && (fRadius >= min) // Wegen Verteilung deren Wertebereich nicht 0.0 bis 1.0 ist.
				   && (fRadius <= max) // Wegen Verteilung deren Wertebereich nicht 0.0 bis 1.0 ist.
				  )	// ende Bedingung
				 bRadiusOkay = true;
				
				
			}
			else
			{
				fRadius = pDistribution->GetRandomValue(); // Zufallswert rauspicken
				fRadius = fRadius * (max - min) + min;
				if(	  (fRadius > 0.0) // Kugel mit keinem oder negativem Radius, ist keine Kugel! Auf wenn minRadius was anderes sagt.
				   && (fRadius >= min) // Wegen Verteilung deren Wertebereich nicht 0.0 bis 1.0 ist.
				   && (fRadius <= max) // Wegen Verteilung deren Wertebereich nicht 0.0 bis 1.0 ist.
				  )	// ende Bedingung
				 bRadiusOkay = true;
				
			}	

		} while( (bRadiusOkay == false) && (nTry++ < nMaxTrys)); // nicht ewig mit dieser Verteilung versuchen
		if(bRadiusOkay == false) // Suche erfolglos
		{	// Letzten Radius als negativen Wert zurueckgeben.
			// Dient nur fuer moegliche Analyse funktionen.
			if(bVerbose)	cout << "Fail with radius = " << fRadius << "; min = " << min << "; max = " << max << endl;
			if(fRadius > 0.0)
				fRadius = 0.0 - fRadius;
		}
	}
	if(bVerbose)	cout << "Radius = " << fRadius << endl;
	return fRadius;
} // END OF double SpongeBob::getRandomRadius(const int nMaxTrys) const
////////////////////////////////////////////////////////////////////////////////




/* ------------------------------------------------------------------------------------ */
void SpongeBob::PeriodicCopy( const sphere s_in, sphere * s_out ){
	/*
	 * creates a set of (in total) 27 spheres forming a periodic arrangement
	 * with respect to the unit cube centered at ( 0., 0., 0. )
	 * 
	 * INPUT:
	 * one sphere
	 *
	 * OUTPUT:
	 * s_out[0] ... s_out[26] (ATTENTION: s_out has to be sized appropriately beforehand!)
	 * 
	 */

// 	const double p_shift[27][3] =
// 				{
// 					{ -0.5,  -0.5,  -0.5},		{  0.0,  -0.5,  -0.5},		{  0.5,  -0.5,  -0.5},
// 					{ -0.5,   0.0,  -0.5},		{  0.0,   0.0,  -0.5},		{  0.5,   0.0,  -0.5},
// 					{ -0.5,   0.5,  -0.5},		{  0.0,   0.5,  -0.5},		{  0.5,   0.5,  -0.5},
// 					{ -0.5,  -0.5,   0.0},		{  0.0,  -0.5,   0.0},		{  0.5,  -0.5,   0.0},
// 					{ -0.5,   0.0,   0.0},		{  0.0,   0.0,   0.0},		{  0.5,   0.0,   0.0},
// 					{ -0.5,   0.5,   0.0},		{  0.0,   0.5,   0.0},		{  0.5,   0.5,   0.0},
// 					{ -0.5,  -0.5,   0.5},		{  0.0,  -0.5,   0.5},		{  0.5,  -0.5,   0.5},
// 					{ -0.5,   0.0,   0.5},		{  0.0,   0.0,   0.5},		{  0.5,   0.0,   0.5},
// 					{ -0.5,   0.5,   0.5},		{  0.0,   0.5,   0.5},		{  0.5,   0.5,   0.5}
// 				};
// 	for(int ip = 0; ip <27; ip++)
// 	{
// 		s_out[ip] = s_in;
// 		s_out.Displace( p_shift[ip] );
// 	}
// }

	//--- sphere newSphere = s_in
	int nSphereOutIndex = 0;
	for(int x = -1; x <= 1; x++)
	{   
		for(int y = -1; y <= 1; y++)
		{	
			for(int z = -1; z <= 1; z++)
			{	
				s_out[nSphereOutIndex].X[0] = s_in.X[0] + (x * fEdgeLength); 
				s_out[nSphereOutIndex].X[1] = s_in.X[1] + (y * fEdgeLength);
				s_out[nSphereOutIndex].X[2] = s_in.X[2] + (z * fEdgeLength);
				s_out[nSphereOutIndex].R = s_in.R; 
				s_out[nSphereOutIndex].coreR = s_in.coreR; 
				nSphereOutIndex++;
			} // next z
		} // next y
	} // next x
} // END OF void SpongeBob::PeriodicCopy( const sphere s_in, sphere * s_out )


/* ------------------------------------------------------------------------------------ */


void SpongeBob::drawVoxelImage(VoxelImage vx_img, unsigned char color, int nSphereCount) //
{	/*if 
	( (nSphereCount <= 0)  && (nSphereCount > MAX_SPHERES) )
	{
		nSphereCount = S_bin;
	}
	for(int i = 0; i < nSphereCount; i++)
	{
		pSpheres[i].DrawVoxelImage(vx_img, color);
	}*/
} // END OF void Spongebob::drawVoxelImage(VoxelImage vx_img, unsigned char col) //


/* ------------------------------------------------------------------------------------ */
